Parameters
filtered_data <- read.csv("data/data-2023-09-11 (2).csv", header = TRUE)
Warning: cannot open file 'data/data-2023-09-11 (2).csv': No such file or directoryError in file(file, "rt") : cannot open the connection
data filtering convertion
mydata_genind
/// GENIND OBJECT /////////
// 698 individuals; 6 loci; 48 alleles; size: 195 Kb
// Basic content
@tab: 698 x 48 matrix of allele counts
@loc.n.all: number of alleles per locus (range: 5-12)
@loc.fac: locus factor for the 48 columns of @tab
@all.names: list of allele names for each locus
@ploidy: ploidy of each individual (range: 2-2)
@type: codom
@call: df2genind(X = filtered_data[, 6:11], sep = "/", ncode = 6, ind.names = filtered_data$indv,
pop = filtered_data$Population, NA.char = "0/0", ploidy = 2,
type = "codom", strata = NULL, hierarchy = NULL)
// Optional content
@pop: population of each individual (group size range: 24-166)
Run basic.stats and render the result
library("pegas")
library(dplyr)
library(tibble)
result <- basic.stats(mydata_hierfstat)
df_resutl_basic<-as.data.frame(result$perloc)
# Weir and Cockrham estimates of Fstatistics - FIS and FST
result_f_stats <- Fst(as.loci(mydata_genind))
result_f_stats <- result_f_stats[,2:3]
colnames(result_f_stats) <- c("Fst (W&C)", "Fis (W&C)")
result_f_stats <- merge(result_f_stats, df_resutl_basic, by="row.names",all.x=TRUE)
colnames(result_f_stats)[10] <- "Fst (Nei)"
colnames(result_f_stats)[12] <- "Fis (Nei)"
result_f_stats <- result_f_stats %>% column_to_rownames(., var = 'Row.names')
result_f_stats_selec <- result_f_stats %>% select(all_of(selected_stats))
result_f_stats_selec
G-statistic
library(hrbrthemes)
# compute the Gstats
result_f_stats <- result_f_stats %>% mutate(GST = 1-Hs/Ht)
result_f_stats <- result_f_stats %>% mutate("GST''" = (n_pop*(Ht-Hs))/((n_pop*Hs-Ht)*(1-Hs)))
result_f_stats
# Plot with linear trend
p2 <- ggplot(result_f_stats, aes(x=GST, y=Hs)) +
geom_point() +
geom_smooth(method=lm , color="red", se=FALSE) +
theme_ipsum()
p2

HW - Panmixia
library("pegas")
hw.test(as.loci(mydata_genind), B = 1000)
Linkage Disequilibrium
LD(as.loci(mydata_genind$tab), locus = c(1, 2), details = TRUE)
a<-LDscan(as.loci(mydata_hierfstat[,2:7]))
LDmap(a)
head(mydata_hierfstat[,2:7])
library("poppr")
pair.ia(mydata_genind, sample = 9)
shuffle df
#shuffled_matrices <- replicate(n_rep, mat[sample(nrow(mat)), ], simplify = FALSE)
shuffled_matrices <- replicate(n_rep, mat[sample(length(mat), replace = FALSE)], simplify = FALSE)
##################
# shuffle only the genotype and add the pop column later for each matrices.
#in the loop?
###############
# Create a list to store the wc
fst_df <- numeric(sequence_length)
fis_df <- numeric(sequence_length)
# Calculate the average for each shuffled matrix
# Iterate through the shuffled matrices
for (i in 1:n_rep) {
# Calculate the statistics for the i-th matrix
#HERE THE COLUMN POP
merged_df <- cbind(level_pop, shuffled_matrices[[i]])
result_f_stats <- wc(shuffled_matrices[[i]])
result_f_stats <- as.data.frame(result_f_stats$per.loc)
# Extract FST and FIS values
fst_values <- result_f_stats$FST
fis_values <- result_f_stats$FIS
print( fst_values)
# Assign values to the data frames
fst_df <- cbind(fst_df, result_f_stats$FST)
fis_df <- cbind(fis_df, result_f_stats$FIS)
}
# Set row names as in result_f_stats
rownames(fst_df) <- rownames(fis_df) <- rownames(result_f_stats)
result_FST <- fst_df[, -1]
fis_df <-fis_df[, -1]
vec <- seq(1, n_rep)
colnames(result_FST) <- colnames(fis_df) <- vec
result_FST[1,]
count (result_f_stats[,1][1] > result_FST[1,] )
count <- sum(result_f_stats[,1][1] > result_FST[1, ])
# Initialize an empty data frame to store the counts
count_df <- data.frame(
Greater = numeric(length(result_FST)),
Smaller = numeric(length(result_FST))
)
# Compare the values in result_f_stats[1] to result_FST for each column
for (col in colnames(result_FST)) {
greater_count <- sum(result_f_stats[1] > result_FST[, col])
smaller_count <- sum(result_f_stats[1] < result_FST[, col])
count_df$Greater[col] <- greater_count
count_df$Smaller[col] <- smaller_count
}
# Print the count data frame
print(count_df)
######################## ########################
LS0tCnRpdGxlOiAiUiBOb3RlYm9vayIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKIyBQYXJhbWV0ZXJzIAoKYGBge3J9CmZpbHRlcmVkX2RhdGEgPC0gcmVhZC5jc3YoImRhdGEtMjAyMy0wOS0xMSAoMikuY3N2IiwgaGVhZGVyID0gVFJVRSkKCnNlbGVjdGVkX3N0YXRzIDwtIGMoIkhvIiwgIkhzIiwgIkh0IiwgIkZpcyAoVyZDKSIsICJGc3QgKFcmQykiLCAiRmlzIChOZWkpIiwgIkZzdCAoTmVpKSIKKQoKbl9yZXA9MTAKCm5fcG9wID0gOAoKc2VxdWVuY2VfbGVuZ3RoIDwtIGxlbmd0aCg2OjExKSAKYGBgCgoKCiMgZGF0YSBmaWx0ZXJpbmcgY29udmVydGlvbgoKYGBge3J9CmxpYnJhcnkoaGllcmZzdGF0KQoKZmlsdGVyZWRfZGF0YSA8LSBkYXRhLmZyYW1lKGluZHYgPSBwYXN0ZShzdWJzdHIoZmlsdGVyZWRfZGF0YSRQb3B1bGF0aW9uLDEsMyksIHJvdy5uYW1lcyhmaWx0ZXJlZF9kYXRhKSwgc2VwPSIuIiksIGZpbHRlcmVkX2RhdGEpCiMgQ3JlYXRlIG15ZGF0YV9nZW5pbmQKcG9wdWxhdGlvbiA8LSBmaWx0ZXJlZF9kYXRhJFBvcHVsYXRpb24KbXlkYXRhX2dlbmluZCA8LSBkZjJnZW5pbmQoCiAgWCA9IGZpbHRlcmVkX2RhdGFbLDY6MTFdLAogIHNlcCA9ICIvIiwKICBuY29kZSA9IDYsCiAgaW5kLm5hbWVzID0gZmlsdGVyZWRfZGF0YSRpbmR2LAogIHBvcCA9IGZpbHRlcmVkX2RhdGEkUG9wdWxhdGlvbiwKICBOQS5jaGFyID0gIjAvMCIsCiAgcGxvaWR5ID0gMiwKICB0eXBlID0gImNvZG9tIiwKICBzdHJhdGEgPSBOVUxMLAogIGhpZXJhcmNoeSA9IE5VTEwKKQpteWRhdGFfZ2VuaW5kCm15ZGF0YV9oaWVyZnN0YXQgPC0gZ2VuaW5kMmhpZXJmc3RhdChteWRhdGFfZ2VuaW5kKQpgYGAKCgoKIyBSdW4gYmFzaWMuc3RhdHMgYW5kIHJlbmRlciB0aGUgcmVzdWx0CgpgYGB7cn0KbGlicmFyeSgicGVnYXMiKQpsaWJyYXJ5KGRwbHlyKQpsaWJyYXJ5KHRpYmJsZSkKCnJlc3VsdCA8LSBiYXNpYy5zdGF0cyhteWRhdGFfaGllcmZzdGF0KQpkZl9yZXN1dGxfYmFzaWM8LWFzLmRhdGEuZnJhbWUocmVzdWx0JHBlcmxvYykKCiMgV2VpciBhbmQgQ29ja3JoYW0gZXN0aW1hdGVzIG9mIEZzdGF0aXN0aWNzIC0gRklTIGFuZCBGU1QgCnJlc3VsdF9mX3N0YXRzIDwtIEZzdChhcy5sb2NpKG15ZGF0YV9nZW5pbmQpKQpyZXN1bHRfZl9zdGF0cyA8LSByZXN1bHRfZl9zdGF0c1ssMjozXQpjb2xuYW1lcyhyZXN1bHRfZl9zdGF0cykgPC0gYygiRnN0IChXJkMpIiwgIkZpcyAoVyZDKSIpCnJlc3VsdF9mX3N0YXRzIDwtIG1lcmdlKHJlc3VsdF9mX3N0YXRzLCBkZl9yZXN1dGxfYmFzaWMsIGJ5PSJyb3cubmFtZXMiLGFsbC54PVRSVUUpCmNvbG5hbWVzKHJlc3VsdF9mX3N0YXRzKVsxMF0gPC0gIkZzdCAoTmVpKSIKY29sbmFtZXMocmVzdWx0X2Zfc3RhdHMpWzEyXSA8LSAiRmlzIChOZWkpIgpyZXN1bHRfZl9zdGF0cyA8LSByZXN1bHRfZl9zdGF0cyAlPiUgY29sdW1uX3RvX3Jvd25hbWVzKC4sIHZhciA9ICdSb3cubmFtZXMnKQpyZXN1bHRfZl9zdGF0c19zZWxlYyA8LSByZXN1bHRfZl9zdGF0cyAlPiUgc2VsZWN0KGFsbF9vZihzZWxlY3RlZF9zdGF0cykpCnJlc3VsdF9mX3N0YXRzX3NlbGVjCmBgYAoKIyBHLXN0YXRpc3RpYwoKYGBge3J9CmxpYnJhcnkoaHJicnRoZW1lcykKCgojIGNvbXB1dGUgdGhlIEdzdGF0cwpyZXN1bHRfZl9zdGF0cyA8LSByZXN1bHRfZl9zdGF0cyAlPiUgbXV0YXRlKEdTVCA9IDEtSHMvSHQpCnJlc3VsdF9mX3N0YXRzIDwtIHJlc3VsdF9mX3N0YXRzICU+JSBtdXRhdGUoIkdTVCcnIiA9IChuX3BvcCooSHQtSHMpKS8oKG5fcG9wKkhzLUh0KSooMS1IcykpKQpyZXN1bHRfZl9zdGF0cwoKIyBQbG90IHdpdGggbGluZWFyIHRyZW5kCnAyIDwtIGdncGxvdChyZXN1bHRfZl9zdGF0cywgYWVzKHg9R1NULCB5PUhzKSkgKwogIGdlb21fcG9pbnQoKSArCiAgZ2VvbV9zbW9vdGgobWV0aG9kPWxtICwgY29sb3I9InJlZCIsIHNlPUZBTFNFKSArCiAgdGhlbWVfaXBzdW0oKQpwMgoKYGBgCgojIE1pc3NpbmcgZGF0YQoKYGBge3J9CiMgTGlicmFyaWVzCmxpYnJhcnkoInBvcHByIikKbGlicmFyeSgiaGVhdG1hcGx5IikKCm1pc3NpbmdfZGF0YSA8LSBpbmZvX3RhYmxlKG15ZGF0YV9nZW5pbmQsIHBsb3QgPSBGQUxTRSkKCiMgTWF0cml4IGZvcm1hdAptYXQgPC0gYXMubWF0cml4KG1pc3NpbmdfZGF0YSkKIyBoZWF0bWFwCnAgPC0gaGVhdG1hcGx5KG1hdCwgCiAgICAgICAgICAgICAgIGRlbmRyb2dyYW0gPSAibm9uZSIsCiAgICAgICAgICAgICAgIHhsYWIgPSAiIiwgeWxhYiA9ICIiLCAKICAgICAgICAgICAgICAgbWFpbiA9ICIiLAogICAgICAgICAgICAgICBzY2FsZSA9ICJjb2x1bW4iLAogICAgICAgICAgICAgICBtYXJnaW5zID0gYyg2MCwxMDAsNDAsMjApLAogICAgICAgICAgICAgICBncmlkX2NvbG9yID0gIndoaXRlIiwKICAgICAgICAgICAgICAgZ3JpZF93aWR0aCA9IDAuMDAwMDEsCiAgICAgICAgICAgICAgIHRpdGxlWCA9IEZBTFNFLAogICAgICAgICAgICAgICBoaWRlX2NvbG9yYmFyID0gVFJVRSwKICAgICAgICAgICAgICAgYnJhbmNoZXNfbHdkID0gMC4xLAogICAgICAgICAgICAgICBsYWJlbF9uYW1lcyA9IGMoIlBvcHVsYXRpb24iLCAiTWFya2VyIiwgIlZhbHVlIiksCiAgICAgICAgICAgICAgIGZvbnRzaXplX3JvdyA9IDgsIGZvbnRzaXplX2NvbCA9IDgsCiAgICAgICAgICAgICAgIGxhYkNvbCA9IGNvbG5hbWVzKG1hdCksCiAgICAgICAgICAgICAgIGxhYlJvdyA9IHJvd25hbWVzKG1hdCksCiAgICAgICAgICAgICAgIGhlYXRtYXBfbGF5ZXJzID0gdGhlbWUoYXhpcy5saW5lPWVsZW1lbnRfYmxhbmsoKSkKKQpwCmBgYAoKCgojIEhXIC0gUGFubWl4aWEKCmBgYHtyfQoKCmxpYnJhcnkoInBlZ2FzIikKCmh3LnRlc3QoYXMubG9jaShteWRhdGFfZ2VuaW5kKSwgQiA9IDEwMDApCmBgYAoKCgojIExpbmthZ2UgRGlzZXF1aWxpYnJpdW0KCmBgYHtyfQoKTEQoYXMubG9jaShteWRhdGFfZ2VuaW5kJHRhYiksIGxvY3VzID0gYygxLCAyKSwgZGV0YWlscyA9IFRSVUUpCgphPC1MRHNjYW4oYXMubG9jaShteWRhdGFfaGllcmZzdGF0WywyOjddKSkKCkxEbWFwKGEpCgoKaGVhZChteWRhdGFfaGllcmZzdGF0WywyOjddKQoKCmxpYnJhcnkoInBvcHByIikKcGFpci5pYShteWRhdGFfZ2VuaW5kLCBzYW1wbGUgPSA5KQoKYGBgCgoKIyBzaHVmZmxlIGRmCgpgYGB7cn0KCiNzaHVmZmxlZF9tYXRyaWNlcyA8LSByZXBsaWNhdGUobl9yZXAsIG1hdFtzYW1wbGUobnJvdyhtYXQpKSwgXSwgc2ltcGxpZnkgPSBGQUxTRSkKc2h1ZmZsZWRfbWF0cmljZXMgPC0gcmVwbGljYXRlKG5fcmVwLCBtYXRbc2FtcGxlKGxlbmd0aChtYXQpLCByZXBsYWNlID0gRkFMU0UpXSwgc2ltcGxpZnkgPSBGQUxTRSkKIyMjIyMjIyMjIyMjIyMjIyMjCiMgc2h1ZmZsZSBvbmx5IHRoZSBnZW5vdHlwZSBhbmQgYWRkIHRoZSBwb3AgY29sdW1uIGxhdGVyIGZvciBlYWNoIG1hdHJpY2VzLiAKI2luIHRoZSBsb29wPyAKIyMjIyMjIyMjIyMjIyMjCgoKIyBDcmVhdGUgYSBsaXN0IHRvIHN0b3JlIHRoZSB3Ywpmc3RfZGYgPC0gbnVtZXJpYyhzZXF1ZW5jZV9sZW5ndGgpCmZpc19kZiA8LSBudW1lcmljKHNlcXVlbmNlX2xlbmd0aCkKCiMgQ2FsY3VsYXRlIHRoZSBhdmVyYWdlIGZvciBlYWNoIHNodWZmbGVkIG1hdHJpeAoKIyBJdGVyYXRlIHRocm91Z2ggdGhlIHNodWZmbGVkIG1hdHJpY2VzCmZvciAoaSBpbiAxOm5fcmVwKSB7CiAgIyBDYWxjdWxhdGUgdGhlIHN0YXRpc3RpY3MgZm9yIHRoZSBpLXRoIG1hdHJpeAogICNIRVJFIFRIRSBDT0xVTU4gUE9QCiAgbWVyZ2VkX2RmIDwtIGNiaW5kKGxldmVsX3BvcCwgc2h1ZmZsZWRfbWF0cmljZXNbW2ldXSkKICByZXN1bHRfZl9zdGF0cyA8LSB3YyhzaHVmZmxlZF9tYXRyaWNlc1tbaV1dKSAKICByZXN1bHRfZl9zdGF0cyA8LSBhcy5kYXRhLmZyYW1lKHJlc3VsdF9mX3N0YXRzJHBlci5sb2MpCiAgIyBFeHRyYWN0IEZTVCBhbmQgRklTIHZhbHVlcwogIGZzdF92YWx1ZXMgPC0gcmVzdWx0X2Zfc3RhdHMkRlNUCiAgZmlzX3ZhbHVlcyA8LSByZXN1bHRfZl9zdGF0cyRGSVMKICBwcmludCggZnN0X3ZhbHVlcykKICAjIEFzc2lnbiB2YWx1ZXMgdG8gdGhlIGRhdGEgZnJhbWVzCiAgZnN0X2RmIDwtIGNiaW5kKGZzdF9kZiwgcmVzdWx0X2Zfc3RhdHMkRlNUKQogIGZpc19kZiA8LSBjYmluZChmaXNfZGYsIHJlc3VsdF9mX3N0YXRzJEZJUykKfQoKIyBTZXQgcm93IG5hbWVzIGFzIGluIHJlc3VsdF9mX3N0YXRzCgpyb3duYW1lcyhmc3RfZGYpIDwtIHJvd25hbWVzKGZpc19kZikgPC0gcm93bmFtZXMocmVzdWx0X2Zfc3RhdHMpCnJlc3VsdF9GU1QgPC0gZnN0X2RmWywgLTFdCmZpc19kZiA8LWZpc19kZlssIC0xXQp2ZWMgPC0gc2VxKDEsIG5fcmVwKQpjb2xuYW1lcyhyZXN1bHRfRlNUKSA8LSBjb2xuYW1lcyhmaXNfZGYpIDwtIHZlYwoKCnJlc3VsdF9GU1RbMSxdIAoKY291bnQgKHJlc3VsdF9mX3N0YXRzWywxXVsxXSA+IHJlc3VsdF9GU1RbMSxdICkKY291bnQgPC0gc3VtKHJlc3VsdF9mX3N0YXRzWywxXVsxXSA+IHJlc3VsdF9GU1RbMSwgXSkKCgoKIyBJbml0aWFsaXplIGFuIGVtcHR5IGRhdGEgZnJhbWUgdG8gc3RvcmUgdGhlIGNvdW50cwpjb3VudF9kZiA8LSBkYXRhLmZyYW1lKAogIEdyZWF0ZXIgPSBudW1lcmljKGxlbmd0aChyZXN1bHRfRlNUKSksCiAgU21hbGxlciA9IG51bWVyaWMobGVuZ3RoKHJlc3VsdF9GU1QpKQopCgojIENvbXBhcmUgdGhlIHZhbHVlcyBpbiByZXN1bHRfZl9zdGF0c1sxXSB0byByZXN1bHRfRlNUIGZvciBlYWNoIGNvbHVtbgpmb3IgKGNvbCBpbiBjb2xuYW1lcyhyZXN1bHRfRlNUKSkgewogIGdyZWF0ZXJfY291bnQgPC0gc3VtKHJlc3VsdF9mX3N0YXRzWzFdID4gcmVzdWx0X0ZTVFssIGNvbF0pCiAgc21hbGxlcl9jb3VudCA8LSBzdW0ocmVzdWx0X2Zfc3RhdHNbMV0gPCByZXN1bHRfRlNUWywgY29sXSkKICBjb3VudF9kZiRHcmVhdGVyW2NvbF0gPC0gZ3JlYXRlcl9jb3VudAogIGNvdW50X2RmJFNtYWxsZXJbY29sXSA8LSBzbWFsbGVyX2NvdW50Cn0KCiMgUHJpbnQgdGhlIGNvdW50IGRhdGEgZnJhbWUKcHJpbnQoY291bnRfZGYpCgojIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMgIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIAoKYGBgCgoKCiMKCmBgYHtyfQoKCmBgYAoKIwoKYGBge3J9CgoKYGBgCgojCgpgYGB7cn0KCgpgYGAKCiMKCmBgYHtyfQoKCmBgYAoKIwoKYGBge3J9CgoKYGBgCgojCgpgYGB7cn0KCgpgYGAKCiMKCmBgYHtyfQoKCmBgYAoKIwoKYGBge3J9CgoKYGBgCg==